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The transition metal-oxygen bond appears prominently throughout chemistry and 
solid-state physics. Many materials, from biomolecules to ferroelectrics to the com- 
ponents of supernova remnants contain this bond in some form. Many of these 
materials' properties strongly depend on fine details of the TM-0 bond and intricate 
correlation effects, which make accurate calculations of their properties very challeng- 
ing. We present quantum Monte Carlo, an explicitly correlated class of methods, 
to improve the accuracy of electronic structure calculations over more traditional 
methods like density functional theory. We find that unlike s-p type bonding, the 
amount of hybridization of the d-p bond in TM-0 materials is strongly dependant 
on electronic correlation. 

I. INTRODUCTION 

Transition metal chemistry is a particularly exciting area of research, with applications 
from astrophysics to biology to potential inexpensive high-efficiency solar cells and high- 
temperature superconductivity. Because of the partially filled d-shell, transition metals 
can form many types of bonds and can also exhibit ferroelectric and ferromagnetic ordering. 
Transition metal oxides are particularly interesting because they are one of the most common 
transition metal complexes, and exhibit most of the above effects. This rich physics is quite 
difficult to describe theoretically, however, since electronic correlation is very strong in these 
materials. Current approximate density functional theories tend to perform quite poorly on 
transition metals, particularly in comparison to its quite good accuracy on elements with s 
and p type bonding. Problematic quantities are not hard to find; they include the dipole 
moment in molecules, binding (or cohesive energies), the lattice constants of perovskites, 
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high pressure behavior, and band gaps/excitation energies. 

Rather than attempting to improve the approximate density functional, quantum Monte 
Carlo (QMC) approaches take a different direction-explicitly treating the electronic corre- 
lation in a wave function based approach, while maintaining reasonable scaling with system 
size. It can be made to scale from 0(1) to 0(N 3 ) in the number of electrons jl|], depending 
on the quantity of interest. QMC attains very low upper-bound energies on medium-sized 
electronic problems (up to thousands of electrons at the time of writing), and has been 
used as a benchmark method on s-p systems 3] . Since it treats the electronic correlation 
explicitly in the many-body wave function, it is a promising method for strongly correlated 
TMO systems. 

The goal of this review is to summarize the current state of the art of QMC as applied 
to TMO's. This is a fairly new field, with few calculations. Most of these calculations have 
benchmarked the method to determine the accuracy that one should expect. This accuracy 
has generally been quite high on most of the quantities studied, particularly for energetics. 
In the course of this work, it has also been determined what trial function (starting guess, 
as explained in the methods section) is necessary to obtain this accuracy. The upper-bound 
property of diffusion Monte Carlo has been critical in this success. By this, we have also 
learned that the electronic correlation in transition metal oxides is entangled with the d-p 
orbital hybridization in these materials. 



II. QUANTUM MONTE CARLO 

The most common flavors of Quantum Monte Carlo that have been used on TMO's are 
variational, diffusion, and reptation Monte Carlo (VMC, DMC, and RMC, respectively). 
We will summarize them here; one can find a more complete review in Ref Another 
flavor, auxiliary field Monte Carlo has been used in a few calculations, but will not be 
discussed here. 

VMC is a direct application of the variational theorem. We write the many-body wave 
function as a function of many-body coordinates R = [r 1; r 2 , . . . , r^J and a set of variational 
parameters P. One then approximates the ground state wave function by minimizing the 
energy expectation value 

E(P) = J ^*(R,P)#^(R,P)dR, (1) 
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assuming that the wave function is normalized. For a complicated variational ansatz such 
as we will introduce later, this integral cannot be evaluated analytically. One can, however, 
evaluate it using Monte Carlo by rearranging the integral to read 

E(P)=/|v&(R,P)| 2 ^ p ^ R. (2) 

Since |^(R, P)| 2 is a probability distribution function, one can sample it using Markov chain 
Monte Carlo and evaluate the energy expectation value as an average over the local energy 
El(R) = . The lowest-energy approximate wave function is then found by minimizing 

the energy. In practice, a combination of energy and the variance of the local energy^] or 
variance only 6j is optimized. 

Many wave functions can be used with VMC, since the only requirement is that one can 
evaluate the wave function and its derivatives quickly. For the work covered in this article, 
we start with a Slater determinant of one-particle orbitals, D, or a linear combination of 
Slater determinants. We then multiply D by the explicitly correlated inhomogeneous Jastrow 
correlation factor e u to obtain the Slater- Jastrow variational wave function De u . We write 

U = J2 u ( r Ui r jii r ij) ( 3 ) 

ijl 

where the lower case indices stand for electronic coordinates, and the upper case indices 
are ionic coordinates. There is considerable choice on how to expand u; for concreteness, 
we show one expansion that performs well enough and has been applied to TMO's. The 
correlation factor is expanded in the Schmidt-Moskowitz form 7J]: 



u{r ih r jI: rtj) = c k a k(ru) + Y c ™ b k( r ij) 

k m 

klm 

where the a k and b k functions are written as 

1 -z{r/r cut ) 

l + P k z(r/r cut y 1 ] 

with different (3k for the different types of functions. The polynomial z(x) = x 2 (6 — 8x + 3x 2 ) 
is chosen so the functions go smoothly to zero at r cut =7.5 bohr. The /3k s and all the 
expansion coefficients c et ,c ee , and c ee% are optimized. If there are multiple determinants, 
their coefficients can also be optimized. We then use the VMC wave function as a trial 
function for RMC or DMC. 
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DMC and RMC are based on the so-called imaginary time Schrodinger equation 

-W£H =iH -E )HK,T), (5) 

which has a steady-state solution <3> , the lowest energy eigenfunction with eigenvalue E 
as long as \I/(R, 0) has a non-zero overlap with $o- All non-steady-state solutions converge 
exponentially to the eigenstate <3>o as T goes to infinity. Transforming to an integral equation, 
we have 

$o(Ri) = Jim J G(Ri, R , r)^ r (R )dR , (6) 

where G is the Green's function of the imaginary time Schrodinger equation and \I/t(Ro) is 
the trial wave function that we obtain from VMC Solving for the exact G for large r is as 
difficult as solving for $ , so we choose some constant small value of r for which we know 
G accurately (for example, see Refs 0, la]), and compound the operations (suppressing the 
r dependence of G): 

$o(R) = hm j G(R, R n ) . . . G(R U Ro)^r(Ro)dR dRi • • • dR n . (7) 

Each application of G is interpreted as a stochastic process, in the same way that the 
diffusion equation can be mapped onto Brownian particles and vice versa (in fact, for a free 
particle, the Hamiltonian is — |V 2 and the simulation is a diffusion process). 

DMC performs a simulation of these random particles for large n. All implementations 
of DMC use a particularly clever importance sampling transformation by multiplying the 
imaginary time Schrodinger equation (Eqn [5]) by the trial function \I/<r(R) and working 
with the time-dependent function \1/t(R)\E'(R, r). Since the time dependence is the same, 
it eventually obtains samples distributed according to the probability distribution function 
Pr oo (R) = $o(R)^t(R)- This transformation improves the efficiency of the calculation by 
several orders of magnitude [3( by using information that we already have about the ground 
state in the form of a trial function. The final probability distribution function can be used 
to evaluate the ground-state energy as follows: 

(E ) = J dR* T (R)$ (R)^^, (8) 

since $o is an eigenstate of H and H can operate forwards or backwards. Any operators 
that do not commute with the Hamiltonian will have expectation values that are biased, 
only becoming unbiased in the limit of = $o- 
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We can remove the error in these operators by using reptation Monte Carlo {9!. [lo|. where 
the random walk is performed in the space of paths: s = [R , Ri, . . . , R n -i, RJ- We sample 
the path probability distribution 

U(s) = ^ T (Ro)G(R , Ri) . . . C(R n _i, R n )^(Rn) (9) 

This can be interpreted in several different ways. If we examine the distribution at Rq, we 
can view the samples of Green's functions as acting on ^(Rn), and therefore Pr (Rq) = 
\1/t(Ro)$o(Ro)- This is the same distribution as we obtain in DMC as the path length 
goes to infinity. Alternatively, since G is symmetric on exchange of the two R coordinates, 
the probability distribution of R n is the same. Finally, we can split the path in two, one 
projecting on ^y(R ), and the other projecting on \l/ T (R n ). We then have 

Pru/2 (Rn/2) = (G(R n / 2 , R n /2-i) • • • C(Ri, R ) x I / t(Ro)) 
x (G(R n / 2) Rn/2+1) • • • C(R n _i, R n )\l/T(R n )) 

= $o(R-/2) 

for n — > 00, which allows us to obtain correct expectation values of operators that do not 
commute with the Hamiltonian. 



III. GEOMETRY OPTIMIZATION 

In TMO materials, it is particularly useful to be able to optimize the geometry of the 
system within QMC The usual way of doing this in mean-field calculations is to calculate 
the forces on the atoms and use one of many minimization routines. Unfortunately, there 
are not yet any reliable methods to calculate the force within diffusion Monte Carlo, despite 
much work in that directionjlO, H, 12, isj]- These methods all require high-accuracy trial 



wave functions, which we usually do not have for transition metals. Thus, with the current 
state of the art, we are only able to optimize a few key degrees of freedom using the total 
energies from DMC calculations and line minimization. Even this must be done carefully 
because of the statistical uncertainty in the DMC energy. What follows is the scheme used 
in the work presented here, which has been found to be quite robust. 

According to Bayes' theorem, given a model M and a set of data D, the probability of 
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the model given the set of data is 

P{D) is an unimportant normalization constant and P(M) is called the prior distribution, 
which we are free to set to reflect the a priori probability distribution on the set of models. 
One usually sets P(M) = 1, the unbiased maximum entropy/least knowledge condition. In 
the case of normally distributed data on a set of points {x±,X2, ...,xn}, 

P(D\M) oc exp[- J2( M ( X *) ~ D(xi)) 2 /2a 2 (xi)}, (11) 

i 

where u{x) is the statistical uncertainty of D(x). 

For example, in the case of bond lengths, we can limit our space of models to M(x) = 
c\ + C2X + C3X 2 , for x close to the minimum bond length. This is equivalent to setting 
the prior distribution equal to one for all quadratic functions and to zero for non-quadratic 
functions. One then calculates several data points D(x) with statistical uncertainties a(x). 
The probability distribution function of the bond length b is then obtained by calculating 
the marginal distribution 

(h) = JS(-c 2 /2c 3 - b)P(D\M)P(M)d Cl dc 2 dc 3 
P[ } J P(D\M)P(M)dc 1 dc 2 dc 3 ' 1 j 

This integral is only three-dimensional, and as such could be calculated by a grid method, 
but it is convenient to calculate it by Monte Carlo, by sampling P(D\M)P(M) and binning 
the bond length. The probability distribution function for the bond length is typically a 
Gaussian function to high accuracy, so it can be described as a mean value with a statistical 
uncertainty. 

To make this scheme more efficient, we would like to calculate QMC energies as far 
away from the minimum as possible while still maintaining accuracy. This is because the 
energy changes much more quickly far from the minimum, which mitigates the stochastic 
uncertainties. That is, the energy scale is larger far from the minimum, so less precision is 
necessary. Thus, we should use a fitting function that is valid as far from the minimum as 
possible, while containing as few parameters as possible. For minimum energy geometries, 
it has been found Q, [l5] that the Vignet or modified Morse potentials are quite good for 
this purpose. 
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IV. APPROXIMATIONS 
A. Pseudopotentials 

In QMC, we can increase the efficiency significantly by using pseudopotentials to replace 
the core electrons with an effective potential. This has the effect of removing the large 
fluctuations near the core, which do not contribute much to the valence electrons' correlation, 
which is the important for chemical properties. This introduces two approximations in the 
technique: first, the pseudopotential itself, and second, the small localization error[]i| in 
diffusion Monte Carlo. 

It has been found that small-core pseudopotentials are necessary for high accuracy on 
transition metals [17I fis|. On the 3d metals, which are the primary focus in this paper, this 
means a Ne-core pseudopotential. The reason for this is that the 3d electrons occupy much 
the same space as the semi core 3p and, to a lesser extent, the 3s electrons. Since the 3d 
electrons are strongly affected by bonding, they in turn interact with the semicore. This 
interaction will change with correlation and chemical environment, so we must include the 
semicore electrons in accurate electronic structure calculations. This is not unique to QMC 



and is generally done in density functional theory where high accuracy is needed 



B. Finite size errors 

When performing calculations for extended systems such as crystals, it is necessary to 
introduce periodic boundary conditions. This is an approximation on two levels. The 
first is the standard one-body level that is corrected by using reciprocal space sampling 
(i.e., k-points). The second level is inherent in a many-body correlated method, where the 
periodic boundary conditions force the electron to interact unphysically with its periodic 
image. This is similar to the finite simulation cell error in classical molecular dynamics 
simulation. This is typically corrected by either modifying the Coulomb interaction to 
remove the spurious interaction [^(J or by an correction 15, 21]. Neither of these methods 
has clearly been demonstrated to be superior, and both methods or similar ones have been 
used successfully. Even with these corrections, a QMC calculation of an extended system 
usually involves on the order of 40 to 100 atoms, regardless of the size of the primitive cell, 
followed by extrapolation to infinite size. 
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C. Fixed node 



The algorithms described above are exact when the wave function can be written as a 
positive function, since then ^t&o is a probability distribution function. For fermions, it is 
not usually the case that has the same zeros as the exact ground state, so we make the 
fixed-node approximation, where the nodal surface of the exact wave function are assumed 
to be the same as the trial wave function. This approximation typically results in recovering 
90-95% of the correlation energy, and can be relaxed, but at the cost of exponential scaling 
of the system sizejsj]. 

Given that the pseudopotential localization approximation is usually quite small for en- 
ergy differences 22], we are mostly concerned with the fixed-node error. The Jastrow factor 
does not change the nodes of the wave function, so in the method outlined above, the nodes 
(and thus the final accuracy) are fixed to be the nodes of the Slater determinant of orbitals 
from the mean-field method. It is currently not feasible to vary the orbital expansion di- 
rectly for a large system, since the number of parameters grows to the thousands for even 
moderately sized systems. However, partial optimizations can be done, and, as we shall see, 
are very effective for transition metal-oxygen systems. 



V. TM-O MOLECULES 



Simple molecular systems are excellent starting points for the study of transition metal 
oxides, since they are small enough to study carefully in a reasonable amount of time, and 
are also treatable by accurate but expensive quantum chemistry techniques like Coupled 
Cluster. This provides an additional much-needed data point to compare accuracy of the 
various electronic structure methods. 



A. Near-optimal one-particle orbitals 



Wagner and Mitas 23j performed the first calculations using DMC on simple two-atom 
transition metal oxides (TiO and MnO), and found a strong dependence of the calculated 
binding energy on the orbitals used in the Slater determinant. They used the B3LYP hybrid 
DFT/Hartree-Fock functional, and varied the percentage of Hartree-Fock mixing. They 
found the optimal percentage to be very close to the semi-empirical value fitted by Becke for 
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his B3PW potential 24] . We have plotted the energy gain of B3LYP orbitals versus Hartree- 
Fock for the first five transition metal monoxide molecules in Fig [TJ Upon examining the 
orbitals, they found a large difference in the d-p hybridization for both TiO (Fig [2]) and 
MnO. This is a direct consequence of the importance of electronic correlation in transition 
metals. 

To understand the importance of the one-particle orbitals, one can conceptually divide 
the total energy in three parts, each described by a different part of the wave function: 

• One-body and antisymmetry: the Slater determinant 

• Two-body electron interaction: Jastrow factor 

• Higher orders : implicit diffusion Monte Carlo wave function. 

The first part, the Slater determinant, determines the nodes of the wave function and there- 
fore the ultimate accuracy of the calculation. Empirically, in materials containing only s and 
p-type elements, these three parts are almost independent of each other-the Hartree-Fock 
orbitals are close to optimal for a Slater- Jastrow wave function. In transition metal oxides, 
however, this situation changes, and the two-body and higher interactions are strong enough 
to change the one-body part significantly. In TMO's, this effect seems to be largely in the 
d-p hybridization between oxygen and the transition metal. 

By using the reptation Monte Carlo algorithm, we can obtain the unbiased one-particle 
density within the fixed- node approximation (Fig [3]), which gives further insight into the 
importance of correlation in the one-particle density. QMC tends to enhance the density 
in the bonding region (the hybridization) over both Hartree-Fock and B3LYP, but is not 
able to completely repair the erroneous Hartree-Fock density because of the fixed-node 
approximation. This is the reason for the large energy gain from using B3LYP orbitals to 
fix the nodal surface. 



B. Energetic Performance 

The total energy of a system is quite important for determination of lowest-energy spin 
states, competing phases, reactions, etc, and is a place where traditional density functional 
theory has encountered difficulties on transition metal oxides. In Table [H we compare the 




ScO TiO VO CrO MnO 

FIG. 1: The energy gain in DMC from using B3LYP orbitals as a function of the metal monoxide. 
The line is a guide to the eye. Taken from Ref [251 ] . 

binding energy obtained by DMC using B3LYP orbitals and several other methods. We find 
excellent accuracy, with the RMS deviations of DMC within the experimental uncertainty for 
most materials. CrO is the only molecule with a large deviation from experiment; however, 
it is not very far outside the experimental uncertainty. DMC is also able to consistently 
obtain a minimum energy bond length with errors below 0.01 A (Table HT|) . better than any 
other published result. 



C. Dipole moments 

While energetics are very important for electronic structure calculations, one is also often 
interested in non-energetic properties, such as dipole moments. There has been little work 
done on such things within QMC, even in the context of simpler s and p systems. To our 



knowledge, the only study of dipole moments other than on TMO's is of the CO molecule 321 ] . 
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FIG. 2: The d-p hybridization orbital (doubly occupied) isosurface for TiO in Hartree-Fock (left) 
and B3LYP (right). B3LYP enhances the hybridization significantly, which leads to lower energy 



in QMC. Figure generated using VMD and POV-Ray 



26, 



27] 




FIG. 3: The density of the Ti-0 molecule projected onto the Ti-0 axis in the bonding region for 
various methods. 
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Method 


ScO 


TiO 


VO 


CrO 


MnO 


RMS 


T A Tool 

LDA 28] 


9.09 


9.13 


O AO 

8.48 


6.26 


6.51 


2.19 


CCSD(T)I29I 


6.71 


6.64 


6.13 


4.20 


3.43 


0.31 


TPSSh[28] 


7.11 


7.18 


6.44 


4.45 


4.62 


0.38 


DMC[25] 


7.06(3) 


6.81(3) 


6.54(3) 


3.98(2) 


3.66(3) 


0.21 


AFQMC[30] 




7.02(21) 






3.79(34) 




Exp [31] 


7.01(12) 


6.92(10) 6.44(20) 4.41(30) 


3.83(8) 






TABLE I: Binding energies of the first five transition metal monoxides by different theoretical 
methods, along with RMS deviations from the experiment (all in eV). Statistical uncertainties in 
units of 10 -2 eV are shown in parentheses for Monte Carlo and experimental results. Zero point 
energy corrections are estimated to be much less than the uncertainty in experiment. There are 
too few AFQMC data to calculate meaningful RMS values. 



Method 


ScO 


TiO 


VO 


CrO 


MnO 


RMS 


LEA [28] 


1.644 


1.597 


1.564 


1.584 


1.602 


0.033 


CCSD(T)[29] 


1.680 


1.628 


1.602 


1.634 


1.66 


0.011 


TPSSh[28] 


1.659 


1.613 


1.582 


1.612 


1.628 


0.012 


DMC[25J 


1.679(2) 1.612(3) 1.587(3) 1.617(4) 1.652(4) 0.008 


Exp[31] 


1.668 


1.623 


1.591 


1.621 


1.648 






TABLE II: Bond lengths in Afor the first five transition metal monoxide molecules. 

A primary reason for this lack of calculations is that until the development of RMC, there 
has not been an easy to implement method to obtain expectation values without the mixed- 
estimator bias. The commonly used methods, pure diffusion Monte Carlo and forward- 



walking [33|, |34J, |35| do not scale well with the system size[36|], since they suffer from increased 
fluctuations of weights as the number of particles increases. One can also use extrapolated 
estimation, where the expectation value of an operator is estimated as (O) = 2((D)dmc ~ 
(O)vmc, but that method introduces an additional approximation that one would like to 
avoid if possible. 

RMC, on the other hand, scales quite well, and is easily applicable to medium-sized 
systems such as TMO molecules. As we have noticed above, the electronic correlation 
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Method 


ScO 


TiO 


VO 


CrO MnO 


T T "\ A Tool 

LDA 28] 


3.57 


3.23 


3.10 


3.41 


CCSD(T)[29j 


3.91 


3.52 


3.60 


3.89 4.99 


TPSSh[28] 


3.48 


3.43 


3.58 


3.97 


RMC [25] 


4.61(5) 


4.11(5) 


4.64(5) 4.76(4) 5.3(1) 


Exp [37] 


4.55 


3.34(1) [38] 


3.355 


3.88 



TABLE III: Dipole moments in Debye. The fixed-node RMC results have been obtained with a 
single determinant of B3LYP orbitals. See text for an analysis of the errors involved for the case 
of TiO. 

and hybridization are very intertwined, and therefore, the electronic correlation and dipole 
moment are also closely related. In Table 1111} we report the dipole moments for the first five 
transition metal monoxides using RMC with B3LYP orbitals. RMC obtains dipole moments 
much higher than that found in experiment, which is somewhat surprising given the high 
accuracy seen in energetic properties. We will explore the fixed node approximation and its 
effect on the dipole moment in the next section. 



D. Beyond the Slater- Jastrow form 

In this section, we explore one of the biggest advantages of the QMC method-the ability 
to go beyond a Slater- Jastrow trial function if needed. As we saw in the previous section, 
RMC with the Slater- Jastrow trial function does not obtain dipole moments in agreement 
with experiment. The dipole moment is very sensitive to electronic correlation, and we wish 
to perform as accurate a calculation as possible to approach the exact value. We can do this 
in QMC by expanding the wave function in determinants. We write the trial wave function 
as 

* T (R)= (Eqa)^, (13) 

where the ZVs are determinants of one-particle orbitals, e u is the Jastrow factor, and q's 
are variational parameters. These determinants and the initial coefficients are taken from a 
Configuration Interaction calculation, and the coefficients are reoptimized using Variational 
Monte Carlo in the presence of the Jastrow factor. 
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This last reoptimization step is crucial, since the DMC energy increases if the CI coeffi- 
cients are kept constant. This is a result of the strong correlation of these systems-the first 
order correlations are taken care of by the Jastrow factor, which the CI calculation tries to 
describe (inefficiently) with determinants. 

In Fig HI we see the convergence of this expansion for TiO. The energy has a smooth 
convergence in the number of determinants, but the dipole moment oscillates significantly, 
with smaller oscillations as the number of determinants increases. The final result is ap- 
proximately 3.8(1) Debye, a significant change from the Slater- Jastrow trial wave function, 
but still quite far from the experimental value of 3.34(1). While this calculation is probably 
not at the exact limit, the dipole moment does not appear to change enough to reconcile 
with experiment. Somewhat reassuringly, though, the Coupled Cluster value also predicts 
a larger value for the dipole moment, so it is possible that the experiment may be in er- 
ror. More studies of non-energy properties using Quantum Monte Carlo are sorely needed, 
however, to obtain an estimate of the expected accuracy. 



VI. SOLIDS 



Calculations on extended TMO systems using QMC are particularly challenging, since 
QMC suffers not only from one body finite size effects (i.e., that described by k-point sam- 
pling), but also from many-body finite size effects, which require large simulation cells. For 
this reason, complete studies as those reported above for molecules are not usually attain- 
able, and most work is still in progress. We will discuss a few preliminary studies and 
a few private communications of work that remains unpublished at the time of this writ- 
ing. Clearly, the details of the calculations may change, so this section is meant more as a 
comment on the current state of the art. 



Using QMC, there have been studies of the antiferromagnet NiO 39|, |40] and MnO 17]. 



Except for Tanaka 



391 ]. who performed a very rough optimization of the lattice constant 



within Variational Monte Carlo, all the published studies calculated only the cohesive energy, 
which comes quite close to experiment (Table HVj) for the materials available. In the very 
recent work of Kolorenc and Mitasj^l], they obtain similar accuracy for the cohesive energy of 
FeO and also obtain the correct ordering of phases for that material, which DFT mispredicts. 
In most of these materials, researchers have found a large dependence on the mean-field 
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FIG. 4: The number of determinants versus the energy and dipole moment for TiO.The dipole 
moments are shifted downwards by 0.1 Debye to correct for the pseudopotential error. 

Material DMC binding energy (eV) Experimental mean-field orbitals 
NiO[40j 9.442(2) 9.5 Hartree-Fock 

MnO[17] 9.40(5) 9.5 

BaTi0 3 [42j 31.2(3) 31.57 LDA 

FeO[4JJ 9.47(4) 9.7 PBE0[43j 

TABLE IV: Cohesive energies for several materials using QMC, all calculated per formula cell. 
Also listed are the optimal mean-field orbitals if reported. LDA is the local density approximation 
of DFT, and PBE0 is a hybrid functional. 

orbitals used, with the optimal orbitals ranging from Hartree-Fock to LDA. Apparently, 
there is no universal optimal mean-field method. 

Wagner and Mitas 42j have also reported using the Bayesian optimization scheme to find 
the minimum energy lattice constant of BaTiC>3, which is well-known to be underestimated 
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by over 1% in the local density approximation to density functional theory, and overestimated 
by a similar margin in the gradient corrections. This 1% error in the lattice constant can 
affect the calculated spontaneous polarization up to 50%, so even this small error is not 
acceptable for a truely first-principles description of this material. DMC obtains a cubic 
lattice constant in error only by 0.015 ± 0.005 A, which is somewhat less than half a percent, 
a significant improvement over the density functional results. Also, in BaTiOa, there is an 
energy gain in DMC of ~1 eV/formula cell by using LDA orbitals instead of Hartree-Fock 
orbitals, and they report that it is due to a similar change in d-p hybridization that is seen 
in the transition metal monoxide molecules. 



VII. CONCLUSIONS 

On the systems that have been tested thus far, QMC offers unprecedented accuracy in a 
completely first-principles and scalable method, particularly in the energetics of the systems. 
The d-p hybridization of transition metal oxides is strongly affected by electronic correlation. 
Using QMC methods, we can clearly see this, both by investigating the minimum-energy 
orbitals and by examining the one-particle density and dipole moment within QMC. The 
dipole moment in particular is strongly affected by the level of correlation present in the 
quantum mechanics approximation. 

On TMO molecules, we have a significant gain in the total energy on expansion into 
determinants, of about 0.5 eV. This means that we are relying on cancellation of errors 
for the high accuracy of QMC, although to a much lesser degree than post-Hartree-Fock 
approaches and DFT. We see this error in the dipole moment, which does not benefit from 
cancellation of errors. On the molecules, however, we can use a brute-force approach by 
expanding in determinants and come quite close to the true ground state. However, this 
kind of expansion will ultimately fail for large systems, since the number of determinants 
grows very quickly with system size. In order to reliably check the QMC results, it is vital to 
develop new reasonably scaling wave functions that go beyond the Slater- Jastrow form. Some 
work has been done in this direction with the RVB44], Pfaffianj^. and backflow 46, 47 ] 



wave functions in QMC. These wave functions' accuracy should be tested on TMO systems in 
the future. Equally important are optimization schemes within VMC that can systematically 
minimize the energy with respect to the wave functions' parameters despite the stochastic 



17 

nature of VMC, which is under serious investigation j^, |4{|. Finally, we need to be able 
to calculate forces within QMC accurately and efficiently The current state of the art is 
not sufficient to treat transition metal oxides[l^], and the Bayesian method of geometry 
optimization is only efficient for a few dimensions. 

The future looks promising for QMC calculations of TMO solids, with the only drawback 
that the calculations are very expensive on today's computers, since one must use a large 
supercell. However, the scaling with system size is quite favorable, and QMC is very easy to 
operate in parallel, so it can take advantage of low-cost processors. It has already been shown 
for a few important transition metal oxide solids that QMC can obtain binding energies and 
other energetic properties with excellent accuracy, well worth the additional cost when high 
accuracy is needed. It remains to be seen how well the method performs for non-energetic 
properties, and what sort of trial wavefunctions are necessary to obtain even higher accuracy. 
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